delta = 0.12;
k = 0;
eigenvalue = 10;

x = sym('x');
y = sym('y');
z = sym('z');
A = [10 0 0; 0 1 0; 0 0 8/3];
R0 = [0; 0; 8*28/3];

u0 = sym('u0');
u1 = sym('u1');
u2 = sym('u2');
u3 = sym('u3');
u4 = sym('u4');
u5 = sym('u5');
u6 = sym('u6');
u7 = sym('u7');
u8 = sym('u8');
u9 = sym('u9');
u10 = sym('u10');

f0 = sym('f0');
f1 = sym('f1');
f2 = sym('f2');
f3 = sym('f3');
f4 = sym('f4');
f5 = sym('f5');
f6 = sym('f6');
%f7 = sym('f7');
%f8 = sym('f8');
%f9 = sym('f9');
%f10 = sym('f10');


u0 = [x - y; y; z];
u1 = -A * u0 - R0     - R1(u0)     -     R2(u0, u0);
u2 = -A * u1 - R1(u1) - R2(u1, u0) -     R2(u0, u1);
u3 = -A * u2 - R1(u2) - R2(u2, u0) - 2 * R2(u1, u1) -      R2(u0, u2);
u4 = -A * u3 - R1(u3) - R2(u3, u0) - 3 * R2(u2, u1) - 3  * R2(u1, u2) -      R2(u0, u3);
u5 = -A * u4 - R1(u4) - R2(u4, u0) - 4 * R2(u3, u1) - 6  * R2(u2, u2) - 4 *  R2(u1, u3) -       R2(u0, u4);
u6 = -A * u5 - R1(u5) - R2(u5, u0) - 5 * R2(u4, u1) - 10 * R2(u3, u2) - 10 * R2(u2, u3) - 5   * R2(u1, u4) -       R2(u0, u5);
%u7 = -A * u6 - R1(u6) - R2(u6, u0) - 6 * R2(u5, u1) - 15 * R2(u4, u2) - 20 * R2(u3, u3) - 15  * R2(u2, u4) - 6   * R2(u1, u5)   - R2(u0, u6);
%u8 = -A * u7 - R1(u7) - R2(u7, u0) - 7 * R2(u6, u1) - 21 * R2(u5, u2) - 35 * R2(u4, u3) - 35  * R2(u3, u4) - 21  * R2(u2, u5)  - 7 * R2(u1, u6)  - R2(u0, u7);
%u9 = -A * u8 - R1(u8) - R2(u8, u0) - 8 * R2(u7, u1) - 28 * R2(u6, u2) - 56 * R2(u5, u3) - 70  * R2(u4, u4) - 56  * R2(u3, u5)  - 28 * R2(u2, u6) - 8 * R2(u1, u7)  - R2(u0, u8);
%u10= -A * u9 - R1(u9) - R2(u9, u0) - 9 * R2(u8, u1) - 36 * R2(u7, u2) - 84 * R2(u6, u3) - 126 * R2(u5, u4) - 126 * R2(u4, u5) - 84 * R2(u3, u6) - 36 * R2(u2, u7) - 9 * R2(u1, u8) - R2(u0, u9);

f0 = R0 + R1(u0) + R2(u0, u0);
f1 = R1(u1) + R2(u1, u0) + R2(u0, u1);
f2 = R1(u2) + R2(u2, u0) + 2 * R2(u1, u1) + R2(u0, u2);
f3 = R1(u3) + R2(u3, u0) + 3 * R2(u2, u1) + 3 * R2(u1, u2)  + R2(u0, u3);
f4 = R1(u4) + R2(u4, u0) + 4 * R2(u3, u1) + 6 * R2(u2, u2)  + 4 * R2(u1, u3)  + R2(u0, u4);
f5 = R1(u5) + R2(u5, u0) + 5 * R2(u4, u1) + 10 * R2(u3, u2) + 10 * R2(u2, u3) + 5 * R2(u1, u4)   + R2(u0, u5);
f6 = R1(u6) + R2(u6, u0) + 6 * R2(u5, u1) + 15 * R2(u4, u2) + 20 * R2(u3, u3) + 15 * R2(u2, u4)  + 6 * R2(u1, u5) + R2(u0, u6);
%f7 = R1(u7) + R2(u7, u0) + 7 * R2(u6, u1) + 21 * R2(u5, u2) + 35 * R2(u4, u3) + 35 * R2(u3, u4)  + 21 * R2(u2, u5) + 7 * R2(u1, u6) + R2(u0, u7);
%f8 = R1(u8) + R2(u8, u0) + 8 * R2(u7, u1) + 28 * R2(u6, u2) + 56 * R2(u5, u3) + 70 * R2(u4, u4)  + 56 * R2(u3, u5) + 28 * R2(u2, u6) + 8 * R2(u1, u7) + R2(u0, u8);
%f9 = R1(u9) + R2(u9, u0) + 9 * R2(u8, u1) + 36 * R2(u7, u2) + 84 * R2(u6, u3) + 126 * R2(u5, u4) + 126 * R2(u4, u5) + 84 * R2(u3, u6) + 36 * R2(u2, u7) + 9 * R2(u1, u8) + R2(u0, u9);
%f10= R1(u10) + R2(u10, u0) + 10 * R2(u9, u1) + 45 * R2(u8, u2) + 120 * R2(u7, u3) + 210 * R2(u6, u4) + 252 * R2(u5, u5) + 210 * R2(u4, u6) + 120 * R2(u3, u7) + 45 * R2(u2, u8) + 10 * R2(u1, u9) + R2(u0, u10);

s0 = s_n(delta, 0, k + eigenvalue);
s1 = s_n(delta, 1, k + eigenvalue);
s2 = s_n(delta, 2, k + eigenvalue);
s3 = s_n(delta, 3, k + eigenvalue);
s4 = s_n(delta, 4, k + eigenvalue);
s5 = s_n(delta, 5, k + eigenvalue);
s6 = s_n(delta, 6, k + eigenvalue);
%s7 = s_n(delta, 7, k + eigenvalue);
%s8 = s_n(delta, 8, k + eigenvalue);
%s9 = s_n(delta, 9, k + eigenvalue);
%s10 = s_n(delta, 10, k + eigenvalue);

j0 = s0 * Q(f0 - k * u0);
j1 = s1 * Q(f1 - k * u1);
j2 = s2 * Q(f2 - k * u2);
j3 = s3 * Q(f3 - k * u3);
j4 = s4 * Q(f4 - k * u4);
j5 = s5 * Q(f5 - k * u5);
j6 = s6 * Q(f6 - k * u6);
%j7 = s7 * Q(f7 - k * u7);
%j8 = s8 * Q(f8 - k * u8);
%j9 = s9 * Q(f9 - k * u9);
%j10= s10 * Q(f10 - k * u10);

i0 = j0 + [1 0 0] * u0;
i1 = j0 + j1 + [1 0 0] * u0;
i2 = j0 + j1 + j2 + [1 0 0] * u0;
i3 = j0 + j1 + j2 + j3 + [1 0 0] * u0;
i4 = j0 + j1 + j2 + j3 + j4 + [1 0 0] * u0;
i5 = j0 + j1 + j2 + j3 + j4 + j5 + [1 0 0] * u0;
i6 = j0 + j1 + j2 + j3 + j4 + j5 + j6 + [1 0 0] * u0;
%i7 = j0 + j1 + j2 + j3 + j4 + j5 + j6 + j7 + [1.44 -1 -1] * u0;
%i8 = j0 + j1 + j2 + j3 + j4 + j5 + j6 + j7 + j8 + [1.44 -1 -1] * u0;
%i9 = j0 + j1 + j2 + j3 + j4 + j5 + j6 + j7 + j8 + j9 + [1.44 -1 -1] * u0;
%i10= j0 + j1 + j2 + j3 + j4 + j5 + j6 + j7 + j8 + j9 + j10 + [1.44 -1 -1] * u0;


fid = fopen('c:\\Lorenz\\polynom_i0.m','wt');
i0 = char(i0);
fprintf(fid, 'function answer1 = polynom_i0(x, y, z)\n answer1=');
fprintf(fid,'%s',i0);
fprintf(fid, ';\n end\n');
fclose(fid);

fid = fopen('c:\\Lorenz\\polynom_i1.m','wt');
i1 = char(i1);
fprintf(fid, 'function answer = polynom_i1(x, y, z)\n answer=');
fprintf(fid,'%s',i1);
fprintf(fid, ';\n end\n');
fclose(fid);

fid = fopen('c:\\Lorenz\\polynom_i2.m','wt');
i2 = char(i2);
fprintf(fid, 'function answer = polynom_i2(x, y, z)\n answer=');
fprintf(fid,'%s',i2);
fprintf(fid, ';\n end\n');
fclose(fid);

fid = fopen('c:\\Lorenz\\polynom_i3.m','wt');
i3 = char(i3);
fprintf(fid, 'function answer = polynom_i3(x, y, z)\n answer=');
fprintf(fid,'%s',i3);
fprintf(fid, ';\n end\n');
fclose(fid);

fid = fopen('c:\\Lorenz\\polynom_i4.m','wt');
i4 = char(i4);
fprintf(fid, 'function answer = polynom_i4(x, y, z)\n answer=');
fprintf(fid,'%s',i4);
fprintf(fid, ';\n end\n');
fclose(fid);

fid = fopen('c:\\Lorenz\\polynom_i5.m','wt');
i5 = char(i5);
fprintf(fid, 'function answer = polynom_i5(x, y, z)\n answer=');
fprintf(fid,'%s',i5);
fprintf(fid, ';\n end\n');
fclose(fid);

fid = fopen('c:\\Lorenz\\polynom_i6.m','wt');
i6 = char(i6);
fprintf(fid, 'function answer = polynom_i6(x, y, z)\n answer=');
fprintf(fid,'%s',i6);
fprintf(fid, ';\n end\n');
fclose(fid);